function [ModDuration, YearDuration, PerDuration] = bnddury(Yield, CouponRate, ...
    Settle, Maturity,varargin)
%BNDDURY Duration of Bond Given Yield.
%   BNDDURY calculates the Macaulay and Modified duration of NUMBONDS fixed
%   income securities given yield to maturity. This function will determine the
%   duration for a bond regardless of whether the bond's coupon structure
%   contains short or long first or last coupon periods (i.e. regardless of
%   whether the coupon structure is synched to maturity). This function will
%   also determine the Macaulay and Modified duration for a zero coupon bond.
%
%   [ModDuration, YearDuration, PerDuration] = bnddury(Yield, CouponRate, ...
%          Settle, Maturity)
%
%   [ModDuration, YearDuration, PerDuration] = bnddury(Yield, CouponRate, ...
%          Settle, Maturity, Period, Basis, EndMonthRule, IssueDate, ...
%          FirstCouponDate, LastCouponDate, StartDate, Face)
%
%   [ModDuration, YearDuration, PerDuration] = bnddury(Yield, CouponRate, ...
%          Settle, Maturity,'Param1','Value1',...)
%
%   Optional Inputs: Period, Basis, EndMonthRule, IssueDate, FirstCouponDate,
%                    LastCouponDate, StartDate, Face, LastCouponInterest,
%                    CompoundingFrequency, DiscountBasis
%
%   Note:
%   - All non-scalar or empty matrix input arguments must be either NUMBONDSx1
%     or 1xNUMBONDS conforming vectors.
%   - Fill unspecified entries in input vectors with NaN.
%   - Dates can be serial date numbers or date strings.
%   - Optional inputs can be specified as parameter value pairs.  If
%     LastCouponInterest, CompoundingFrequency or DiscountBasis are input,
%     optional inputs must be specified as parameter value pairs.
%     Otherwise, optional inputs may be specified by order according to the
%     help.
%
%   Inputs:
%        Yield - Yield to maturity on a semi-annual basis.
%
%   CouponRate - Coupon rate in decimal form.
%
%       Settle - Settlement date.
%
%     Maturity - Maturity date.
%
%   Optional Inputs:
%            Period - Number of coupons payments per year.
%                     Possible values include:
%                     0, 1, 2 (default), 3, 4, 6, 12
%
%
%             Basis - Day-count basis.
%                     Possible values include:
%                     0 - actual/actual (default)
%                     1 - 30/360 SIA
%                     2 - actual/360
%                     3 - actual/365
%                     4 - 30/360 PSA
%                     5 - 30/360 ISDA
%                     6 - 30/360 European
%                     7 - actual/365 Japanese
%                     8 - actual/actual ISMA
%                     9 - actual/360 ISMA
%                    10 - actual/365 ISMA
%                    11 - 30/360 ISMA
%                    12 - actual/365 ISDA
%                    13 - bus/252
%
%      EndMonthRule - End-of-month rule; default is 1 (in effect)
%                     0 - Rule is NOT in effect for the bond(s)
%                     1 - (default) Rule is in effect for the bond(s) (meaning
%                         that a security that pays coupon interest on the last
%                         day of the month will always make payment on the last
%                         day of the month)
%
%         IssueDate - Bond issue date.
%
%   FirstCouponDate - Irregular or normal first coupon date.
%
%    LastCouponDate - Irregular or normal last coupon date.
%
%         StartDate - Forward starting date of payments.
%
%              Face - Face value of the bond; default is 100.
%
%   LastCouponInterest - Compounding convention for computing the yield of
%                        a bond in the last coupon period, i.e.: with only
%                        the last coupon and the face value to be repaid.
%                        Choices are 'simple' or 'compound'.
%
%   CompoundingFrequency - Compounding frequency for yield calculation.  By
%                          default, SIA bases (0-7) and BUS/252 use a semi-annual
%                          compounding convention and ISMA bases (8-12) use
%                          an annual compounding convention.
%
%   DiscountBasis - Basis used to compute the discount factors for
%                   computing the yield.  The default behavior is for SIA
%                   bases to use the actual/actual day count to compute
%                   discount factors, and for ISMA day counts and BUS/252
%                   to use the specified basis.
%
%   Outputs:
%   ModifiedDuration - [NUMBONDSx1 vector] of the Modified duration.
%
%       YearDuration - [NUMBONDSx1 vector] of the Macaulay duration in years.
%
%        PerDuration - [NUMBONDSx1 vector] of the Periodic Macaulay  duration.
%
%   See also BNDDURP, BNDCONVY, BNDCONVP.

%   Copyright 1995-2009 The MathWorks, Inc.
%   $Revision: 1.7.2.13 $   $Date: 2009/11/05 16:58:59 $

% Checking input arguments
if nargin < 4
    error('Finance:bnddury:tooFewinputs', 'Too few inputs.');
end

% Check to see whether we have the case of ordered inputs or PV pairs
if ~isempty(varargin)
    if ischar(varargin{1})
        
        p = inputParser;
        
        p.addParamValue('discountbasis',NaN);
        p.addParamValue('compoundingfrequency',NaN);
        p.addParamValue('lastcouponinterest',{'compound'});
        
        p.addParamValue('period',[]);
        p.addParamValue('basis',[]);
        p.addParamValue('endmonthrule',[]);
        p.addParamValue('issuedate', []);
        p.addParamValue('firstcoupondate', []);
        p.addParamValue('lastcoupondate', []);
        p.addParamValue('startdate', []);
        p.addParamValue('face',[]);
        
        try
            p.parse(varargin{:});
        catch ME
            newME = MException('finance:bnddury:optionalInputError',...
                'Error in input arguments');
            newME = addCause(newME,ME);
            throw(newME)
        end
        
        CompFreq = p.Results.compoundingfrequency;
        Basis =  p.Results.basis;
        DiscountBasis = p.Results.discountbasis;
        
        try
            [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate, Settle, Maturity,[], Basis);
        catch ME
            throwAsCaller(ME)
        end
    else
        CompFreq = NaN;
        DiscountBasis = NaN;
        
        try
            [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate, Settle, Maturity,varargin{:});
        catch ME
            throwAsCaller(ME)
        end
    end
else
    
    CompFreq = NaN;
    DiscountBasis = NaN;
    
    try
        [CouponRate, Settle, Maturity,~,Basis] = instargbond(CouponRate, Settle, Maturity);
    catch ME
        throwAsCaller(ME)
    end
end

% The scalar expansion done inside instargbond may not be correct since it
% doesn't consider "Yield".
try
    [CouponRate, Settle, Maturity, Yield, CompFreq, DiscountBasis] = ...
        finargsz(1, CouponRate, Settle, Maturity, Yield, CompFreq, DiscountBasis);
catch ME
    throwAsCaller(ME)
end

% Call cfamounts for bond parameters
[CFlowAmounts, ~, TFactors] = cfamounts(CouponRate, Settle, Maturity, varargin{:});

%Get the size of the porftolio
nBonds = size(CFlowAmounts, 1);

% Figure out the compounding frequency here
if isnan(CompFreq)
    CompFreq = 2*ones(nBonds,1);
    if isnan(DiscountBasis)
        i = isisma(Basis);
    else
        i = isisma(DiscountBasis);
    end
    CompFreq(i) = 1;
end

%-----------------------------------------------------------------------------
% Assign all intermediate variables
TF = TFactors(:, 2 : end);
CF = CFlowAmounts(:, 2 : end);
YLD = Yield(:, ones(1, size(CF, 2)));
PeriodMat = CompFreq(:, ones(1, size(CF, 2)));

%-----------------------------------------------------------------------------
% Calculate periodic Macaulay duration

A = nan .* ones(nBonds, 1);
B = A;

for i = 1:nBonds
    
    a = TF(i, :).*CF(i, :) ./ ((1 + YLD(i, :) ./ PeriodMat(i, :)).^TF(i, :));
    
    SumInd = ~isnan(a);
    
    A(i, 1) = sum(a(SumInd), 2);
    
    b = CF(i, :) ./ ((1 + YLD(i, :) ./ PeriodMat(i, :)).^TF(i, :));
    
    SumInd = ~isnan(b);
    
    B(i, 1) = sum(b(SumInd), 2);
end

PerDuration = A ./ B;

% Calculate yearly duration
YearDuration = PerDuration ./ CompFreq;

% Now calculate modified duration
ModDuration = YearDuration ./ (1 + Yield./CompFreq);